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I review, in a personal perspective, the history of the theory of non- 
relativistic bound states in QED and QCD from the Bethe-Salpeter 
equation to the construction of effective field theories. 



1.1. Introduction 

The study of bound states and, in particular, of non-rclativistic bound 
states has accompanied the quantum theory from its beginning through all 
its subsequent turning points up to what is now the Standard Model of 
particle physics. At the beginning it was the description of the hydrogen 
atom that led to the foundation of quantum mechanics, later the Lamb 
shift contributed to the development of relativistic field theories and renor- 
malization, which eventually led to the foundation of Quantum Electro- 
dynamics (QED); similarly, in the seventies, quarkonium played a special 
role in the foundation of Quantum Chromodynamics (QCD). The special 
role of non-relativistic bound states in particle physics is due to the strik- 
ing experimental signatures that they provide and the fact that analytical 
(perturbative) methods are able to describe the relevant features of these 
signatures. 

Despite this, it has proven very difficult to carry out theoretical analyses 
of a precision comparable with the data, in part due to the high quality 
of the data, but largely owing to the difficulties in performing bound-state 
calculations. These may be traced back to the presence of different energy 
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scales that make it a challenge to maintain a consistent book-keeping in the 
calculations. 

Let us consider a non-relativistic particle of mass m that propagates 
in a potential V (in the case of a Coulomb potential: V = —a/r). If the 
momentum of the particle is non relativistic, then p ^ mv, u <C 1 being 
the velocity of the particle. In the threshold region, the velocity is such 
that mv'^ ^ V. The balance between kinetic energy and potential creates 
the bound state: the particle propagator G cannot be computed order by 
order in V, but comes from resumming all potential insertions in the free 
propagator Gq ~ 1/{E — p'^/2m): 



The function G = !/(£' — /2m — V) exhibits poles in correspondence of 
the bound-state energies En ^ mv^ (= —ma^/i-n? in the Coulombic case, 
which implies u ~ a and l/r ~ ma), the residues at the poles, 0* satisfy 
the equation: 



which is the Schrodinger equation for a non-relativistic bound state whose 
wave function is 0„. 

Hence, the non-relativistic dynamics of a particle close to threshold is 
characterized by a hierarchy of energy scales: m 3> mv 3> mv^ . The scale 
of the mass is sometimes called "hard" , the scale of the typical momentum 
transfer, or inverse size of the system, mu, is called "soft" and the scale 
mv^ is called "ultrasoft" . 

At the level of non-relativistic quantum mechanics, m does not play any 
dynamical role, for the kinetic energy is p^ /2m rather than -y/p^ + m^. The 
contributions from the other scales are accounted for by the Schrodinger 
equations (jl.ip or (|1.2p . The solutions of the Schrodinger equation are 
non-relativistic bound states of typical energy of order mv"^ and typical 
momentum (or inverse size) of order mv. 

One may expect that a more complicated picture will emerge in a rela- 
tivistic field theory, although the leading dynamics should be still described 
by a Schrodinger equation. In a relativistic field theory description of the 
bound state, we will have, besides the bound state, other degrees of free- 
dom, for instance photons (in QED) and gluons (in QCD) emitted and 
exchanged by the bound state; for each of them, modes associated to each 
of the energy scales, to, mv and mv"^ will appear. We shall discuss bound 
states in relativistic field theories in the next section. 



G — G*o ^" G(^VG. 



(1.1) 




(1.2) 
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1.2. The Bethe Salpeter equation 

Let us consider a particle and an antiparticle (e.g. an electron and a 
positron or a quark and an antiquark) that interact near threshold. In the 
centre-of-mass frame, their momenta p and energies E are small compared 
to their masses m: p/m ^ v 1. We assume that we may express the 
interaction perturbatively in terms of Feynman diagrams. This is always 
the case in QED, but does not need to be so in QCD where, at the typical 
hadronic scale AqcD) perturbation theory breaks down. Non-relativistic 
bound states in QCD are made by heavy quarks: this means that at least 
the quark mass is larger than Aqcd • A bound state of a heavy quark and a 
heavy antiquark is called quarkonium (examples are charmonium, a charm- 
anticharm bound state, and bottomonium, a bottom-antibottom bound 
state; a top-antitop bound state, which would be toponium, has no time to 
form due to the rapid top quark weak decay, however, near threshold, the 
bound-state enhancement should be visible in the top-antitop production 
cross section). A perturbative treatment of quarkonium, which requires 
mw, mv^ ^ Aqcd, is justified only for top-antitop pairs near threshold and 
possibly for the ground state of bottomonium. 



1 

E — p^ jm — V 



as (1 -I- as/v + . . .) 

Fig. 1.1. Rcsummcd propagator near threshold. 

How does the bound state emerge in a near threshold interaction? For 

certain sets of graphs, like those in Fig. 11.11 the perturbative expansion 

breaks down when as v (for definiteness, we will consider here and in 

the following figures the QCD case: continuous lines stand for quarks and 

antiquarks, and the curly lines for gluons; the strong coupling constant is 

as). The summation of all as/v contributions leads to the appearance of a 

bound-state pole of order mv'^ ~ ma^ in the resummed propagator. Indeed, 

in the leading non-relativistic limit, when the quark/antiquark propagators 

i 1 ± 7° 

can be approximated by - — — — ^-7- — and the gluon 

^ ±p° + E/2-p^/2m + ie 2 ^ 

exchange by — (close to threshold we may expand in |g^|/|q| ^ 1; 7^ is a 




April 6, 2009 12:52 



World Scientific Review Volume - 9in x 6in 




4 



Antonio Vairo 



Dirac matrix) the Green's function shown in Fig. [TTT] satisfies Eq. (|f .ip . 

Beyond the leading non-relativistic Umit, diagrams wiU be much more 
compUcated to calculate and contributions from the different energy scales 
will get entangled. This happens for any diagram, but the annihilation 
diagram shown in Fig. 11.21 provides a rather immediate way to see it. As- 
suming that the incoming quarks are near threshold, the different gluons 
entering the diagram are characterized by different scales: the annihilation 
gluons have a typical energy of order m; binding gluons carry the momen- 
tum of the incoming quarks, which is of order mv, and ultrasoft gluons, 
sensitive to the intermediate bound state, have energies of the order of the 
binding energy, i.e. mv^. 



Fig. 1.2. Annihilation diagram contributing to the quarkonium decay width. 

The entanglement of the different energy modes makes it difficult to 
organize a full relativistic calculation of the bound state in QED or QCD. 
An equation suitable for bound states in field theory, formally similar to 
was suggested almost sixty years ago by H. Bethe and E. Salpeter:^ 



where G is the two-particle Green's function, Gq the product of the free 
propagators of the two particles and the kernel K is the sum of all ampu- 
tated irreducible two-particle diagrams. Equation (jl.3p does not represent 
an expansion, because, like (jl.ip . a bound state emerges only from the sum 
of all interactions, at least those shown in Fig. 11.11 However, unlike (jl.ip . 
the Bethe-Salpeter equation is not homogeneous in the momentum scale 
and an exact solution is unknown. To make the Bethe-Salpeter equation 
useful, the strategy has been to isolate from K a kernel Kc containing the 
leading contribution responsible for the formation of the bound state, i.e. 




mv 



G = Gq + GqKG, 



(1.3) 
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the Coulomb potential, and expand around it (see, for instance, Ref. 2). 
The most refined approach in this strategy can be found in Ref. 3 (see also 
Ref. 4) : Kc is chosen in such a way that the corresponding Bethe-Salpeter 
equation, Gc = Go + Goi^cGc, may be solved in an analytically closed form, 
and the full Green's function expanded around the exact solution: 



where 5K = K — Kc- Since, in the non-relativistic limit, Kc becomes the 
Coulomb potential, Gc is nothing else than a relativistic modification of the 
solution of Eq. p.ip for V equal to the Coulomb potential, which is known 
since long time.^ The difference between Eq. (11.31) and Eq. (|1.4p is that 
the latter is a perturbative expansion in the kernel while the former is not. 

The Bethe-Salpeter equation was the only systematic tool to treat 
bound states in field theory until the end of the eighties. However, around 
that time, it became increasingly clear that perturbative calculations for 
QED bound states, to which the Bethe-Salpeter equation had been mostly 
applied, could not be push beyond the reached limit if not at the cost of 
a formidable amount of work. It shows the difficulty of the approach the 
fact that going from the calculation of the ma^ correction in the hyper- 
fine splitting of the positronium ground state^ to the ma^ In a term^'^ took 
twenty-five years! The main problem was the lack of an efficient way of 
disentangling the contributions coming from the different energy scales and 
organize them in a perturbative expansion (techniques for asymptotic ex- 
pansions of Feynman integrals near threshold would be developed later^): 
each Feynman diagram would contribute to the observables with a series in 
the coupling constant. No obvious counting rules were available even for the 
leading term of the series. Also gauge invariance did not provide a useful 
organizational tool, since it was very cumbersome to isolate gauge-invariant 
subsets of diagrams in K}'^^^^ 

In the late seventies and eighties, systematic calculations of quarkonium 
observables started (for a recent review see Ref. 14). The complicated dy- 
namics of QCD made it more apparent that a treatment based on the 
Bethe-Salpeter equation was inadequate to perform high-precision quarko- 
nium calculations. First, not all of the quarkonium scales are in general per- 
turbative, lower ones may not be, so that a separation of scales is necessary 
to achieve factorization. Second, even if a perturbative treatment would 
be possible (like for the bottomonium ground state and for tt threshold 
production), the number and topology of diagrams makes the calculation 
prohibitive. It was felt that somehow going back to the Schrodinger equa- 



G = Gc + GcSKG, 



(1.4) 
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tion and identifying a quarkonium potential would lead to a more treatable 
problem. In Refs. 15-19, a quarkonium potential was derived from the 
quark-antiquark scattering amplitude. In the same years, focusing in par- 
ticular on toponium and tt threshold physics, a similar program was carried 
out by W. Kummer and coUaborators^""^^ (see also the Ph.D. thesis in Ref. 
25). In this case, the starting point was the Bethe-Salpeter equation and 
the generalization to QCD of the solution of the Bethe-Salpeter equation 
for positronium found in Ref. 3. Still, the goal was not the solution of the 
Bethe-Salpeter equation itself, but the derivation of a potential, facing, in 
the process, some of the problems that, in a few years, would have led to 
(and found a solution with) the construction of effective field theories for 
non-relativistic bound states. Among the problems mentioned or addressed 
at that time were the infrared sensitivity of the potential, the inclusion of a 
finite decay width (in Refs. 22,23, one can find addressed, for the first time 
in a formal way, how to include the top-quark instability beyond leading 
order), gauge invariance. The infrared sensitivity of the potential will be 
discussed in Sec. 11.51 



1.3. NRQED/NRQCD 

In QED and QCD, one may take advantage of the hierarchy of scales that 
characterizes non-relativistic bound states by expanding Green's functions 
in the ratios of low energy scales over large energy scales. Working out 
these expansions, however, turns out to be cumbersome and does not lead 
to a straightforward and easy way to organize the calculation. If such an 
expansion is instead implemented at the Lagrangian level, it leads to the 
construction of an effective field theory (EFT) . In the effective field theory 
the large scale is integrated out from the beginning and does not appear 
anymore in the Green's function. The terms in the EFT Lagrangian are 
organized as an expansion in powers of the inverse of the large scale that 
has been integrated out leading to a straightforward power counting. 

The first EFT introduced for non-relativistic bound states in QED and 
QCD has been non-relativistic QED/QCD (NRQED/NRQCD). The large 
scale that is integrated out in NRQED/NRQCD is the mass m of the bound- 
state constituents. The degrees of freedom of NRQED/NRQCD are non- 
relativistic fermions and antifermions, and photons/giuons of energy and 
momentum smaller than m; they build up the operators 0„ of the La- 
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grangian. The Lagrangian is organized as an expansion in 1/rn: 

^NRQED/NRQCD = ^ C„(as(m), ^) X (1-5) 
n 

Since, once On has been run down to energies lower than m, the expecta- 
tion value of On scales like mv or smaller scales, Eq. (|1.5p provides, for 
any physical observable, a perturbative expansion in the ratio of the scale 
mv or smaller scales over m. The Wilson coefficients c„ are non analytical 
in the scale m and function of the factorization scale fi. They are calcu- 
lated by equating, "matching", amplitudes in QED/QCD with amplitudes 
in NRQED/NRQCD order by order in l/m and in the coupling constant 
since in both theories we have that a, as{m) ^ 1. The matching may be 
performed on scattering amplitudes, hence in a manner completely inde- 
pendent of the bound state. This is not surprising: the formation of the 
bound state takes place at a scale, mv, which is much smaller than m. 




QCD NRQCD 
Fig. 1.3. Matching to NRQCD. 



The diagram in Fig. 11.21 corresponds, via the optical theorem, to the 
imaginary part of the diagram shown on the left in Fig. 11.31 The same 
process would be described in NRQCD by the diagram shown on the right 
in Fig. 11.31 i.e. by a diagram where the two hard gluons coming from the 
annihilation are replaced by a contact interaction. The difference between 
the two diagrams is compensated by the Wilson coefficient c ~ as{m)'^. 

As our example may suggest, NRQCD is particularly well suited to 
describe heavy quarkonium decay and production. ^'''^^ It is in the theory 
of quarkonium production that NRQCD has perhaps achieved its major 
success by explaining, in the nineties, the quarkonium production data at 
the Tevatron by a new mechanism allowed by the symmetries of NRQCD, 
the octet mechanism, but missed by previous approaches (see Ref. 14 and 
references therein). 

Applications of NRQED have started in the nineties and with time have 
led to many new results (for some early works, see Refs. 29-31). However, 
the progress in high precision calculations in NRQED/NRQCD has been 
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slowed down by two major shortcomings: first, the fact that soft and ultra- 
soft degrees of freedom still remain entangled in NRQED/NRQCD, second, 
the use in early NRQED/NRQCD calculations of a cut-off regularization 
scheme. The first difficulty led to a power counting that was non homo- 
geneous and to perturbative calculations that still involved two scales. To 
overcome this difficulty, lower energy EFTs were developed; we shall dis- 
cuss some of them in Sec. 11.51 The second difficulty, on one hand, pushed 
the development of lattice NRQCD"^^ (see Ref. 33 for recent results on the 
bottomonium spectrum), on the other hand, addressed analytical studies 
towards a consistent formulation of NRQCD in dimensional regularization. 

1.4. The bound state in dimensional regularization 

Surprisingly, it was only few years after NRQCD had been introduced that 
an EFT for mesons made of a single heavy quark, the heavy quark effec- 
tive theory (HQET), was formulated. '^^ In the two-fermion sector, the La- 
grangian of HQET contains the same operators as the NRQCD Lagrangian. 
However, HQET is a quite different theory from NRQCD: HQET contains 
only a single dynamical scale, Aqcd, which governs its power counting. As 
a consequence, the kinetic energy, which is of order AQQj~,/m, is suppressed 
with respect to the binding energy, which is of order Aqcd, while, in a 
non-relativistic bound state, the two are of the same order. 

It is precisely because, in HQET, propagators are expanded in the ki- 
netic energy that we may use dimensional regularization in loop calcula- 
tions. This has led to a rapid, vast and very successful use of the HQET in 
precision studies of D and B mesons. ^'"^ Instead, keeping the kinetic energy 
in the denominators of the propagators, as the power counting of NRQCD 
seems to suggest, turns out to be disastrous and leads to the break down 
of the power counting. The reason is that, in dimensional regularization, 
integrals are not cut-off at high momenta and hard scale poles are going to 
contribute if present in the denominators. Once, this had been realized in 
Ref. 36, it became also clear that the way out was to compute the matching 
to NRQCD in the same way as the matching to the HQET, i.e. order by 
order in 1/m. Since both in NRQCD and in the HQET the matching con- 
ditions are computed in the same way, the two Lagrangians are the same: 
not only the operators of the two theories coincide in the two-fermion sec- 
tor, but also their matching coefficients do. Obviously, in order to compute 
observables with the NRQCD Lagrangian, the usual non-relativistic power 
counting rules, different from the HQET ones, should be used. 



April 6, 2009 12:52 World Scientific Review Volume - 9in x 6in KummerMemoriarhep 



9 



Having understood how to treat the bound state in dimensional regular- 
ization, opened, finally, the doors to analytical high-precision calculations 
also for non-relativistic bound states in NRQED/NRQCD. 

1.5. pNRQED/pNRQCD 

The problem of disentangling the soft from the ultrasoft scale in 
NRQED/NRQCD was addressed immediately after dimensional regular- 
ization was established as an useful tool for non-relativistic bound state 
calculations also. The history and details of the developments that have 
ultimately led to the construction of EFTs for the ultrasoft degrees of free- 
dom of NRQED/NRQCD have been recollected in Ref. 37 and we refer the 
interested reader to it. Here, we would like just to stress the importance 
that the process of tt production near threshold (see Ref. 14 and refer- 
ences therein) has played in these developments, providing the only near 
threshold, heavy quark-antiquark system in nature entirely accessible in 
perturbation theory. As it was mentioned before, this very special feature 
of the tt system near threshold had already been appreciated by the groups 
working on the subject at the beginning of the nineties and, in particular, 
by the Vienna group. 

In the following, in order to illustrate some general features, we will con- 
centrate on the EFTs for ultrasoft degrees of freedom of NRQED /NRQCD 
known as potential NRQED^^.ss and potential NRQCD''°'''i (for an aher- 
native formulation see Ref. 42 and the review in Ref. 43). The large scale 
that is integrated out in pNRQED/pNRQCD is the typical momentum 
transfer of the bound state, which, in coordinate space, is associated with 
the inverse of the typical distance r between the two heavy particles. The 
degrees of freedom of pNRQED/pNRQCD are non-relativistic fermions and 
antifermions, and photons/gluons of energy and momentum smaller than 
mv. They build up the operators Ok,n of the Lagrangian; the operators 
may be also chosen to be explicitly gauge invariant. The Lagrangian is 
organized as an expansion in 1/m, inherited from NRQED/NRQCD, and 
an expansion in r (multipole expansion) , which is characteristic of the new 
EFT: 

'^^pNRQED/pNRQCD = X! ^ Cfe(as(m),^) X F„(r,^,^') X r"Ofe,„(^')- 

k,n 

(1.6) 

Since, once Ok,n has been run down to the lowest energy mu^, the expec- 
tation value of Ok^n scales like mw^, Eq. (|1.6[) provides, for any physical 
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observable, a perturbative expansion in the ratio of mv'^ over mv or to. 
The Wilson coefficients Ck are those inherited from NRQED/NRQCD, the 
Wilson coefficients Vn are the new ones of pNRQED/pNRQCD. They are 
non analytical in the scale r and function of the new factorization scale 
yCt'. They are calculated by matching, order by order in r. Green's function 
in NRQED/NRQCD with Green's function in pNRQED/pNRQCD. In pN- 
RQED, the matching may be also done order by order in a. In pNRQCD, 
as{mv) <§; 1 holds only for tightly bound states (short-range quarkonia, 
e.g. the bottomonium ground state or tt near threshold), while, in general, 
higher excited quarkonium states (long-range quarkonia) are not accessible 
by perturbation theory. This means that we can rely on an expansion in 
as{mv) only for the former states, while for the latter states the matching 
has to be done in a non-perturbative fashion. 




+ 



1 



E ~ p'^/m — V{r, fi, ii') 
Fig. 1.4. Matching to pNRQCD. 

Let us consider the NRQCD diagram shown in Fig. II. 31 and the part of it 
where soft gluons are exchanged between the quark and antiquark. The sum 
of all soft-gluon exchanges would be described in pNRQCD by the diagram 
shown on the right in Fig. 11.41 where the single line stands for a quark- 
antiquark propagator in a color-singlet configuration, 1/{E — p^/m — V), 
the double line for a quark-antiquark pair in a color-octet configuration, the 
curly line for ultrasoft gluons and the circle with a cross for a chromoelec- 
tric dipole interaction ~ ^^r • E (f> that comes from multipole expanding the 
gluon fields in the NRQCD Lagrangian. The non-analytical behaviour in r 
of the NRQCD diagram is reproduced in pNRQCD by the Wilson coefficient 
V. Since V, together with p^ /m, makes up the pole of the quark-antiquark 
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propagator, the interpretation of V is obvious: V is the potential describ- 
ing the interaction in the heavy quark-antiquark pair. At leading order in 
the multipole expansion, when we neglect diagrams involving ultrasoft glu- 
ons, the equation of motion of a non-relativistic fermion-antifermion pair 
is nothing else than the Schrodinger equation (|1.2p . 

The Schrodinger equation is the equation governing non-relativistic 
bound states in quantum mechanics. The full rclativistic description pro- 
vided by field theory, which is richer and much more complex, is given by 
the Bethe-Salpeter equation. This complexity arises from the entanglement 
of different energy scales. Once the contributions of all these scales have 
been separated/factorized, we are left with an EFT of the ultrasoft degrees 
of freedom. The Schrodinger equation naturally emerges as the equation of 
motion of these ultrasoft degrees of freedom. But, because the EFT con- 
tains all the richness and complexity of the field theory, although unfolded 
in a systematic and organized way, the Schrodinger equation, which we have 
gotten from the EFT, is much more than the Schrodinger equation of quan- 
tum mechanics we have started with. First, the EFT provides a proper, 
field theoretically founded, definition of the potential: the potential is the 
Wilson coefficient of the dimension six operator of the EFT, containing two 
fermion and two antifermion fields, that encodes all contributions coming 
from modes whose energies and momenta are larger than the binding en- 
ergy. It undergoes renormalization, develops scale dependence and satisfies 
renormalization group equations, which, in perturbation theory, allow to 
resum potentially large logarithms. Moreover, the EFT accounts also for 
effects that cannot be cast in a Schrodinger equation and that are due to the 
coupling of the fermion-antifermion pair with the other ultrasoft degrees of 
freedom. 







1 1 
1 1 
1 1 
1 1 

1 1 
1 1 
1 1 
1 1 







Fig. 1.5. QCD diagrams responsible for the infrared sensitivity of the static potential. 

In QCD, ultrasoft effects affect also the static potential. As first ob- 
served in Ref. 44, they come from the "non-Abelian Lamb shift" ^'^ diagrams 
displayed in Fig. 11.51 At fixed order in perturbation theory, the diagrams 
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are infrared divergent; at order a^, the leading logarithmic correction is' 



The result shows clearly the non-physical nature of the potential, which 
depends on the renormalization scale fi' . The potentially large logarithms, 
ln(r/i'), have been resummed by means of renormalization group equations 
in Ref. 46; subleading corrections have been calculated in Ref. 47. In 
physical observables, like the static energy or the quarkonium mass, the 
scale dependence of Eq. (jl.7p cancels against ultrasoft contributions coming 
from the second diagram in the pNRQCD part of Fig. 11.41 

Higher-order terms in the relativistic expansion may be computed sys- 
tematically in the EFT. Again, the full complexity and symmetries of the 
underlying field theory are not lost in the expansion. So, for instance, 
relativistic invariance imposes specific constraints on the Wilson coeffi- 
cients/potentials of the EFT,'**''*^ which can be tested on the lattice. ^"'^^ 

Applications of pNRQCD and, more in general, of EFTs for the ultrasoft 
degrees of freedom of NRQCD have led to a plethora of new results in 
quarkonium physics (see Refs. 14,37,52-55 for some recent reviews) and, 
in particular, in tt threshold production (see Refs. 56-59 for the present 
status of the art). Also QED calculations have remarkably benefitted from 
the EFT approach and corrections of very high order in perturbation theory 
have been calculated in the last years for many observables after decades of 
very slow or no progress. As an example, we mention that for the hyperfine 
splitting of the positronium ground state the terms of order a®, In^ a and 

In a have been calculated (for recent reviews on positronium precision 
studies and further references we refer to Refs. 60,61). 

In Fig. 11.61 we summarize the hierarchy of EFTs for bound states in 
QED and, for heavy quarks, in QCD. 

1.6. Outlook 

The history of non-rclativistic bound states in the quantum theory had in 
the last century a peculiar spiral behaviour. It started with the Schrodinger 
equation of the hydrogen atom and seemed to have written its ultimate 
chapter with the Bcthe-Salpeter equation in the fifties. However, in face of 
the enormous difficulties in treating bound states in field theory by means 
of the Bethe-Salpeter equation, a long journey started in the seventies 
that took us back to the Schrodinger equation. This coming back, how- 
ever, was not like closing a circle, it was more like building up a spiral. The 
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LONG-RANGE 
QUARKONIUM 



SHORT-RANGE 
QUARKONIUM / QED 



QCD/QED 



m 



perturhative matching 



pertiirbative matching 



mv 



NRQCD/NRQED 



non— perturhative 
matching 



perturhative matching 



pNRQCD/pNRQED 



Fig. 1.6. EFTs for bound states in QED and, for heavy quarks, in QCD. 

Schrodinger equation we have come back to, encompasses all the complexity 
of the Bethe-Salpeter equation, all the richness of field theory, in the elegant 
and systematic setting of non-relativistic effective field theories. The count- 
ing rules and structure of the EFTs have allowed us to perform calculation 
with unprecedented precision, where higher-order perturhative calculations 
were possible, and to systematically factorize short from long range contri- 
butions where observables were sensitive to the non-perturbative, infrared 
dynamics of QCD. 

Non-relativistic EFTs have become nowadays the standard tool to treat 
non-relativistic bound states. Besides QED bound states and quarkonium, 
these include hadronic atoms like pionium,^^ nucleon-nucleon systems, ^'^'^^ 
non-relativistic bound states at finite temperature^^^^^ and many others. 
The modern history of non-relativistic bound states is far from being fi- 
nished and still needs to be told in its full extent. 
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